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Abstract 

In this paper we review the trajectory of a model proposed by 
Stauffer and Weisbuch in 1992 to describe the evolution of the immune 
repertoire and present new results about its dynamical behavior. Ten 
years later this model, which is based on the ideas of the immune net- 
work as proposed by Jerne, has been able to describe a multi-connected 
network and could be used to reproduce immunization and aging ex- 
periments performed with mice. Besides its biological implications, 
the physical aspects of the complex dynamics of this network is very 
interesting per se. The immunization protocol is simulated by intro- 
ducing small and large perturbations (damages), and in this work we 
discuss the role of both. In a very recent paper we studied the aging 
effects by using auto-correlation functions, and the results obtained 
apparently indicated that the small perturbations would be more im- 
portant than the large ones, since their cumulative effects may change 
the attractor of the dynamics. However our new results indicate that 
both types of perturbations are important. It is the cooperative ef- 
fects between both that lead to the complex behavior which allows to 
reproduce experimental results. 
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1 Introduction 



The main task of the immune system is to protect the organism against 
dangerous elements: antigens (virus, bacteria, poison, cell residues, 
etc). Depending on the antigen the immune system may elicit dif- 
ferent kinds of responses: the cell-mediated immune response or the 
humoral response. The models discussed in this paper are related to 
the humoral responses generated by B cells (one of the main classes 
of lymphocytes), which are the cells that produce the antibodies. The 
antibodies produced by a given population of B cells are replicas of its 
molecular receptor. Each molecular receptor exhibited by a given B 
cell population recognize different recognition sites (epitopes) of the 
antigen by lock-key interactions. 

The clonal selection theory is the most accepted theory about 
the B cells and was proposed by Burnet in 1959. It states that the anti- 
gen chooses by pattern recognition the clones of B cells (population 
of B cell and antibodies) that will proliferate. In order to recognize 
any foreign (or dangerous) element the immune repertoire must be 
complete. According to estimates the human immune system is able 
to express at least the order of 10 6 different receptors 1 . Due to 
the completeness of the repertoire, the immune system would be able 
to recognize and be recognized (recognizing epitopes of its own an- 
tibodies), therefore the same mechanism of recognition should work 
for both antibody-antigen and antibody-antibody reactions. In 1974 
Jerne 2 , taking into account these different mechanisms, suggested 
that when the antigen is presented to the organism it will activate a 
set of clones of B cells, leading to the production of specific antibodies; 
those antibodies will in turn be recognized by a second set of clones 
activating them, and so on. Due to the interplay of the mechanisms 
of activation and suppression, this chain of reactions will be finite, 
preventing the "information" from "percolating" through the entire 
system. Therefore the immune system functions as a functional net- 
work 2 , with complexity comparable to the nervous system. Since 
its proposal, only few evidences support the existence of the immune 
network theory E] , and those evidences suggest that if the network 
exists only 20% of the lymphocytes will be activated, while the rest 
of the clones will form a pool of immunocompetent lymphocytes that 
are able to recognize any antigen. 

In what follows we define in section 12 a model introduced by Stauf- 
fer and Weisbuch in 1992 5 , which incorporates the main concepts of 
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Jerne's immune network theory. In section |3] we review the main re- 
sults obtained in the literature, focusing on both the physical jSHZj and 
biological aspects jSj , since among many other attempts to model what 
an immune network could be, this was the first that could reproduce 
the behavior of a real immune system. In section |1] we present new 
results about the short and long time behaviors of the multi-connected 
network under small and large perturbations. We also discuss the re- 
sults obtained about the aging effects in a very recent work [H] on the 
light of these new results. In our concluding remarks (section |SJ) we 
address the future prospects from the physical and biological points 
of view. 

2 The Model 

In 1992 Stauffer and Weisbuch (SW) introduced a cellular automata 
model to describe the evolution of the B cell repertoire. This model 
was inspired in a previous one (based on a difference equation ap- 
proach) introduced by de Boer, Segel and Perelson (BSP) ^0] using 
a shape-space formalism to describe a large-scale immune network. 
Using a discrete shape-space formalism SW associate each point of a 
d-dimensional space to a different molecular receptor (or a clone of B 
cell). In this way each receptor is characterized by d properties, and 
its neighbors in the shape space will correspond to molecular recep- 
tors that differ from it by one of these properties (according to the 
estimates based on a continuous approach if the notion of shape 
space is relevant, then d > 5). 

To each receptor they associate a three-state cellular automaton 
B(f,t) that will describe the concentration of the population charac- 
terized by this receptor f over the time: low (B(f,t) = 0), interme- 
diate (B(f,t) = 1) and high (B(f,t) = 2). Following the ideas of the 
immune network theory, the interactions among different populations 
are based on complementarity (lock-key interaction) and are defined 
by deterministic rules. The time evolution of the cellular automaton 
is based in a deterministic non-local rule: population B(f,t) at site r 
is influenced by the populations at site — r (its mirror image or com- 
plementary shape) and its nearest-neighbors {—f+5r) (representing 
defective lock-key interactions). The influence on the population at 
site r caused by its complementary populations is described by the 
field h(f,t): 
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h{r,t) 



(1) 



r 'E(— r+<5r) 



where for a given r the sum runs over the complementary shape f ' = 
— r and its nearest neighbors. Due to the finite number of states of the 
population B, the maximum value of the field h{f) is h max = 2(2d+l). 
The rules are based on a window of activation for each population 
which is inspired in a log bell shaped proliferation function associated 
to the receptor cross- linking involved in the B cell activation OH] • 
There is a minimum field necessary to activate the proliferation of the 
receptor populations (0±), but for high doses of activation (greater 
than 62) the proliferation is suppressed. The updating rule may be 
summarized as: 



but no change is made if it would lead to B = —1 ox B = 2>. We define 
the densities of sites in state i at time t as Bi(t) (i = 1, 2, 3). 

The initial configurations are randomly generated according to 
the parameter control x, which determines the initial concentrations: 
Bi(0) = -62(0) = x/2, while the remaining L d {\ — x) sites are initiated 
with B(r,0) = 0. 



3 Results obtained from previous stud- 
ies 



The above definition of the model corresponds already to the modified 
version which has been proposed in Ref. SW have shown that for 
the original version of the model [5], there is a stable-chaotic transi- 
tion on the behavior of the discrete system for d > 4, when varying 
the parameter x. As pointed out by the authors, according to the 
characteristics of those behaviors, none of them would be appropriate 
to describe the evolution of the immune repertoire. 

One of us has studied how the system attains the chaotic regime, 
by studying the transient times and periods of the original model close 
to the transition 12 . It was observed that close to the transition, the 
system is trapped in cycles with large but finite periods. In 1995 




if 0i < 
otherwise 




(2) 



4 



Zorzenon dos Santos and Bernardes f| introduced the modified ver- 
sion of SW model described above and have extensively studied the 
behavior of this model in the parameter space. Differently from SW, 
they obtained |H] a stable-chaotic transition for d > 2, and mapped 
all the behaviors on "phase" diagrams for different set of parameters, 
showing that there was a broad transition region between stable and 
chaotic behaviors. 

In the following work 7 u Bernardes and Zorzenon dos Santos have 
investigated the dynamical behavior of the model on this broad tran- 
sition region. They have observed an aggregation-disaggregation dy- 
namics, with clusters splitting and fusing along the time, as a multi- 
connected network. The authors have also studied the behavior of the 
system in the transition region when subjected to antigen presenta- 
tion. The antigen presentation was simulated by flipping randomly 
chosen populations from the inactive state to the highly activated 
state, reproducing in this way the activation of the populations caused 
by the presence of an antigen. The results obtained by the standard 
spread of damage procedure indicate that after the antigen pre- 
sentation, the perturbation first increases and then relaxes after a few 
time steps, indicating that some information about this perturbation 
is incorporated by the system. These results suggested that this model 
would be a good candidate to model a real immune system behavior. 
Zorzenon dos Santos and Bernardes then applied it [S] to reproduce 
immunization and aging experiments performed with mice under mul- 
tiple antigen presentations ^UE3EI] obtained by the group of Nelson 
Vaz. 

In the experiment of immunization Jl] six mice of the same linage 
are subjected to the following protocol: the antigen presentations are 
produced by intra-venal injections of ova. The time interval between 
the first and second presentation is 14 days, but after that the time 
interval between consecutive presentations is fixed to 7 days. Before 
each antigen presentation (except for the first one) they measure the 
amount of specific antibodies in the blood. Among many experiments 
reported by this group, Zorzenon dos Santos and Bernardes chose the 
one showing a refractory (saturation) behavior related to the immu- 
nization protocol Jl] that could not be explained by clonal selection 
theory. However, the aging experiment chosen |15l I16j is one among 
others obtained by the group showing the same result: a reduction on 
the intensity of the response as the system ages. 

In order to simulate the immunization protocol [Hj, for each sam- 
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pie they let the system evolve from its initial configuration to 1000 
time steps, and associate this time step to the birth of the simulated 
mouse. Then random small perturbations, in size and location in the 
shape-space (using different time intervals between them), are pro- 
duced in order to simulate the noise to which the immune system is 
subjected due to the environment, food ingestion, etc. By using the 
arbitrary scale of 1 day corresponding to 5 time steps, it is possible to 
simulate independent antigen presentations without over-exciting the 
system j7]. The immunization protocol is simulated by large pertur- 
bations (one order of magnitude greater than the small ones) being 
produced periodically every 7 days. The protocol was simulated for 
different samples (10) of young and old mice (8 and 24 weeks, respec- 
tively). The results could be interpreted on the light of the immune 
network theory, the refractory behavior being associated to the sat- 
uration of the multi-connected network in incorporating information 
about the perturbation, and the aging effects related to the loss of 
plasticity of the system. 

Since there are only few evidences that support the existence of the 
immune network, this model played a very important role in making 
the connection between some experimental results that could not be 
explained by clonal selection theory, and the immune network theory. 
Moreover the theoretical results indicate that under multiple antigen 
presentation the saturation is due not only to the number of specific 
antibodies produced (as shown in the Elisa experiments) but also the 
ability to incorporate information by including new populations in the 
network. Since the number of populations belonging to the network 
does not change significantly during the time evolution, in order to 
add information, part of the information already present in the net- 
work should be lost. Therefore the profile of expressed antibodies 
would change after each antigen presentation, reflecting such changes. 
This is an aspect that comes out of the theoretical results and should 
be investigated by the experimentalists when the appropriated tools 
become available. The techniques currently available do not allow to 
identify the difference among the specific antibodies. That would be a 
requirement to verify whether they change after multiple presentations 
as suggested by the model. However, there are some experimental ev- 
idences that it may happen |T7j . 

Once the complex multi-connected network generated by the model 
can be used to reproduce the behavior of a real immune system, it 
becomes interesting to investigate its dynamical properties from the 
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physical point of view. This investigation will allow to learn about the 
dynamics of a (real) complex system operating out of equilibrium. The 
characterization of such behavior started with the work of Bernardes 
and Zorzenon dos Santos |18j . who analyzed the distribution of cluster 
sizes and the distribution of the permanence time (i.e. the time inter- 
val during which each population remains activated or belong to the 
multi-connected network) . They found a characteristic cluster size and 
a power law behavior for the distribution of permanence times. The 
characteristic cluster size is associated to the loss of plasticity while 
the power law distribution indicate that the populations belonging to 
the network have no typical permanence time, reflecting the fact the 
dynamical memory is generated by incorporating information about 
the different antigens presented to the system. 

In a very recent work in collaboration with D. Stariolo f£] we have 
used the auto-correlation functions (a common tool to study aging 
effects in e.g. glassy systems) in order to analyze the similarities be- 
tween the aging effects in the multi-connected network and in glassy 
systems. While in glassy systems the aging effects results from the 
loss of plasticity generated by frustration mechanisms, in the network 
these effects are caused by multiple perturbations, since in order to 
incorporate information the system has to adapt to (or satisfy) the 
mechanisms of activation and suppression of the dynamics. The usual 
procedure adopted to study glassy systems consists in taking a "pic- 
ture" of the system at time t w and calculating the number of sites 
that do not change during the following time steps. In this sense 
the auto-correlation function is equal to Ctot(t,t w ) = 1 — hd(t,t w ) = 
1 — HD(t,t w )/N, where hd(t,t w ) [HD(t,t w )] is the normalized (non- 
normalized) Hamming distance between the configurations at times t 
and t w . According to the results obtained, the system without any 
perturbation is driven to a long period attractor after a long transient 
time (10 4 — 10 5 time steps). All the biologically relevant effects are ob- 
served in the transient time of the system. In the purely deterministic 
case, the transient is simply the time it takes for the system to reach 
the attractor (~ 10 4 time steps). When subjected to random small 
perturbations, however, the very notion of a transient becomes fuzzy: 
results in Ref. [2] show that the system approaches a cycle, but is de- 
flected from it by the small perturbations after ~ 10 3 — 10 4 time steps. 
Therefore, the small perturbations will cause the system to change at- 
tractors from time to time due to their cumulative effects. This is 
reflected in the decrease of the auto-correlation functions, as can be 
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Figure 1: Densities vs. time (upper panel) and correlations vs. t — t w (lower 
panel) for a system with small perturbations, L = 50, x = 0.26 and t w = 10 5 . 
C to t corresponds to the autocorrelation function defined in the text. C 2 2 
measures the autocorrelation only in the subspace of sites with B(f, t w ) = 2. 



seen Fig. ^ where the system leaves the cycle it had approached until 
t w = 10 5 time steps (these results are discussed in detail in Ref. [9J. 
The analogue observed in glassy systems corresponds to changes of 
the system to different minima of the energy landscape during the 
relaxation time. Contrary to one's initial intuition that the large per- 
turbations would accelerate the de-correlation process, the changes 
induced by large perturbations do not lead to this effect. Curiously, 
the large perturbations alone (which correspond to the immunizations) 
lead to a much weaker (slower) de-correlation: this is due to the fact 
that they are produced always at the same sites in shape space. Small 
perturbations can be more easily absorbed by the system than the 
large ones since they involve only local changes. The results from 
the autocorrelation functions therefore indicate that the system first 
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Figure 2: Mean Hamming distance vs. time for the complete immunization 
protocol adopted in Ref. jHj and for different ages. Average over 10 runs. 
Standard deviations are not shown for the sake of clarity. 



approaches a cycle (phase space compression), but after some time 
the small perturbations eventually deflect it from the trajectory to its 
"natural" attractor 



4 Results 

This scenario provides a possible explanation for the results regard- 
ing the spread of damage in the system, which managed to reproduce 
the experimental data from immunization experiments with mice (HJ. 
Fig. [2] shows the time evolution of the Hamming distance between a 
system which undergoes only small perturbations, and its copy, which 
undergoes the same small perturbations and the large periodic im- 
munizations. All curves saturate after some time, just like the Elisa 
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Figure 3: Stationary Hamming distance vs. age at which the immunizations 
start. Bars correspond to standard errors (average over 10 runs). 



measurements [Sj. More important, however, is the fact that the value 
of the saturation depends on the age at which immunization started. 
As the age increases, the saturation value decreases, a result which 
was also observed in the mice and is interpreted as a sign of loss of 
plasticity with age [S]. This might be related to the fact that the sys- 
tem is in its transient, while trying to reach its long-period attractor. 
The "older" the system is, the less configurations available there are. 

If we now repeat the computer experiment of immunization (with 
small and large perturbations) for ages greater than ~ 10 3 time steps, 
we obtain a different result: for older systems the stationary Hamming 
distance may increase, as shown in the rightmost curve of Fig. |2J We 
argue that the cumulative effects of the small perturbations could 
be responsible for the change of behavior in the rightmost curve of 
Fig. [2j for sufficiently old systems, the "noise" of the small pertur- 
bations become more important than earlier in the transient, and the 
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Figure 4: Normalized Hamming distance between configurations at times t 
and t + 1, taking into account only activated (B = 2) sites. 



monotonicity of the stationary value of the HD with the age at which 
immunization starts, no longer holds. This result is summarized in 
Fig. |31 It is interesting to point out that this breakdown of the refrac- 
toriness induced by age, takes place on time scales which are compa- 
rable (typically larger than) the life time of the mice. 

The large perturbations, on the other hand, propagate to the entire 
system and some of the active populations generated during the anti- 
gen presentation will remain as part of the network. This mechanism 
depends on the plasticity of the system or, in other words, will be con- 
trolled by the characteristic cluster size of the system |18j . When we 
combine small and large perturbations, the system is driven to a cycle 
that incorporates information about the antigen, but will eventually 
be deviated from that trajectory as soon as the cumulative effects of 
the small perturbations dominate the dynamics In order to ob- 
serve the instantaneous effects of both types of perturbations, we have 
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analyzed the normalized Hamming distance between the same system 
at two consecutive time steps, focusing only on the subspace of sites 
with B(r,t) = 2. This quantity is denoted by HD22(t,t + 1) and re- 
ferred to simply as "instantaneous Hamming distance". It is plotted 
in Fig. HJas a function of time for different situations. 

We start the analysis of Fig. 0] by noticing the behavior for the 
purely deterministic case (lower curve): the plot indicates that typ- 
ically 2.6% of the active network (sites with B{f, t) = 2) changes at 
each time step. The end of the transient is clearly seen at t ~ 15000 
time steps (for this particular realization): the curve then becomes 
periodic, with a very long period. Note that the periodicity of the 
dynamics can be confirmed by checking that the Ct Q t(t w ,t) = 1 for 
some t > t w and for sufficiently large t w (see Fig. ^and Ref. To 
illustrate the periodicity of HD22 (t, t + 1), we plot in Fig. Elthe return 
map of the time series from t = 20000 to t = 25000 time steps (upper 
plot). The extremely long period might lead one to suspect, based 
on the return map, that the system might be chaotic. This suspicion 
is ruled out when one observes the return map from t = 25000 to 
t = 30000 time steps (lower plot). We challenge the reader to find a 
discrepancy between the two return maps. 

Returning to Fig. |11 notice that the size of the changes induced by 
the small perturbations changes near the end of the transient time. 
While the perturbed system no longer attains the cycle, it can still 
feel the effects of the attractor: when close to the cycle, the small per- 
turbations induce larger changes in the active network at each time 
step. If subjected only to the large immunizations, the instantaneous 
Hamming distance take typically larger values. This might suggest 
a contradiction with the results in Ref. !) . according to which the 
system only with large immunizations de-correlates slower. But the 
contradiction is only apparent: HD(t w ,t) measures a long term ef- 
fect, which is dominated by the small perturbations. Fig. |1J on the 
other hand, is a measure of the instantaneous change occurring in the 
system, i.e. a short time effect. Finally, one can combine both types 
of perturbation, the result being even larger instantaneous changes. 

Another aspect of Fig. |1] that should be noted is an interesting 
interplay between the perturbations. In what follows, we refer to the 
typical values of the four curves in Fig. 0] (for large t), from bottom to 
top, as d ("deterministic"), s ("small"), I ("large") and Is ("large and 
small") — therefore d < s < I < Is, as discussed above. The fact that 
Is — s > I — d indicates that the effect of the large immunizations (i.e. 
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HD 22 (t,t+l): Return map for first half of the time series 
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Figure 5: Top panel: return map of HD 2 2(t,t + 1) from t = 20000 to t = 
25000. Bottom panel: return map of HD 22 (t,t + 1) from t = 25001 to 
t = 30000. 



the instantaneous change it induces) is boosted by the presence of the 
small perturbations. The reverse is obviously also true: the previous 
inequality implies Is— I > s—d, indicating that the small perturbations 
induce larger changes in the network if the large immunizations are 
present. This signals a cooperative interaction between two kinds of 
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perturbations which are significantly different in nature, leading to the 
necessary complexity of behavior which allow the model to reproduce 
experiments performed with real immune systems. 

5 Concluding Remarks 

In this paper we have reviewed the main results from the biological 
and physical point of view of a cellular automata model introduced by 
SW in 1992. Since it incorporates the main features of the immune 
network theory and models the functioning of a multi-connected net- 
work, it emerged as a candidate to reproduce the basic dynamics of 
what the real immune network could be. The model was then used 
to model real experiments performed with mice about immunization 
and aging. The immunization experiments are simulated by producing 
small and large perturbations or damages. In this paper we comple- 
ment the results of a very recent work |5] with the study of the role of 
each kind of perturbation. The results obtained in the previous work, 
using auto-correlation functions, apparently suggested that the small 
perturbations would be more important than the large ones from the 
dynamical point of view, since their cumulative effects could change 
the attractor of the dynamics after a few thousand time steps. How- 
ever this study was focusing in the long term behavior of the system. 
Here we have studied its short term behavior and shown that both 
types of perturbation are important: it is their cooperative effects 
that generate the necessary complexity which allowed to reproduce 
experiments performed with real immune systems. 

Despite the understanding of the dynamics we have achieved up to 
this point, there are still some points that should be investigated: what 
happens when changing the frequency of presentation of the small 
perturbations for a fixed value of the large damages, and how these 
results will change when varying the large damage size? This study 
(already in progress) will draw a scenario that would allow us to try to 
apply this model to reproduce other experimental results concerning, 
for instance, cross-reactivity mechanisms observed in immunization 
protocols, and tolerance. 
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